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' Abstract: We perform a comprehensive expforatfon of the Constramed MSSM parameter space 

' employing a Markov Chain Monte Carfo technique and a Bayesian analysis. We compute superpart- 

■ ner masses and other collider observables, as well as a cold dark matter abundance, and compare 
fH , them with experimental data. We include uncertainties arising from theoretical approximations 
Q-i, as well as from residual experimental errors of relevant Standard Model parameters. We delineate 
Oh! probability distributions of the CMSSM parameters, the collider and cosmological observables as 
pH I well as a dark matter direct detection cross section. The 68% probability intervals of the CMSSM 

parameters are: 0.52 TeV < TOi/2 < 1.26 TeV, mo < 2.10 TeV, -0.34 TeV < Aq < 2.41 TeV 

■ and 38.5 < tan/3 < 54.6. Generally, large fractions of high probability ranges of the superpartner 
5_j ■ masses will be probed at the LHC. For example, we find that the probability of mg < 2.7 TeV 

' is 78%, of niq^ < 2.5 TeV is 85% and of m^± < 0.8 TeV is 65%. As regards the other observ- 

ables, for example at 68% probability we find 3.5 x 10~^ < BR{Bs < 1-7 x 10~*, 

1.9 X 10-^° < (5a^uSY < 9 9 ^ iQ-io and 1 x lO^i^^pb < cr|^ < 1 x lO^^pb for direct WIMP detec- 
tion. We highlight a complementarity between LHC and WIMP dark matter searches in exploring 
the CMSSM parameter space. We further expose a number of correlations among the observables, 
in particular between BR{Bs fJ-^fJ-^) and BR{I3 Xgj) or ct^-^. Once SUSY is discovered, this 
and other correlations may prove helpful in distinguishing the CMSSM from other supersymmetric 
models. We investigate the robustness of our results in terms of the assumed ranges of CMSSM 
parameters and the effect of the {g — 2)^ anomaly which shows some tension with the other ob- 
servables. We find that the results for mo, and the observables which strongly depend on it, are 
sensitive to our assumptions, while our conclusions for the other variables are robust. 
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1. Introduction 

Two of the most challenging questions facing particle physics today are the instability of 
the Higgs mass against radiative corrections (known as the "fine-tuning problem") and 
the nature of dark matter (DM). Unlike in the Standard Model (SM), both find plausible 
solutions in the framework of weak scale softly broken supersymmetry (SUSY). Firstly, 
the fine-tuning problem is addressed via the cancellation of quadratic divergences in the 
radiative corrections to the Higgs mass. Secondly, assuming i2-parity, the lightest super- 
symmetric particle (LSP) is a leading weakly interactive massive particle (WIMP) candi- 
date for cold DM (CDM). On the other hand, despite these and other attractive features, 
without a reference to grand unified theories (GUTs), low energy SUSY models suffer from 
the lack of predictivity due to a large number of free parameters {e.g., over 120 in the 
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Minimal Supersymmetric Standard Model (MSSM)), most of which arise from the SUSY 
breaking sector. At present, experimental constraints on superpartner masses from direct 
SUSY searches at LEP and the Tevatron remain fairly mild, although a dramatic improve- 
ment is expected once the LHC comes into operation in 2007. Indirect limits from bounds 
on CP-violation and flavor changing neutral currents are generally much stronger (except 
for the 2nd-3rd generation mixings) but even these can be evaded, for example if SUSY 
breaking is "universal". 

The MSSM with one particularly popular choice of universal boundary conditions at 
the grand unification scale is called the Constrained Minimal Supersymmetric Standard 
Model (CMSSM) Q. The CMSSM is defined in terms of five free parameters: common 
scalar (mo), gaugino {mi/2) and tri-linear (^0) rnass parameters (all specified at the GUT 
scale) plus the ratio tan/3 of Higgs vacuum expectation values and sign(/i), where fi is the 
Higgs/higgsino mass parameter whose square is computed from the conditions of radiative 
electroweak symmetry breaking (EWSB). The economy of parameters in this scheme makes 
it a useful tool for exploring SUSY phenomenology. In addition to experimental limits on 
Higgs and superpartner masses, and strong bounds on SUSY contributions to BR{B — > 
Xs^) and the anomalous magnetic moment of the muon {g — 2)^, a very strong constraint 
limiting the mass parameters of the model from above is provided by the relic abundance 
of the LSP. Within the CMSSM the neutral LSP is the (bino-like) hghtest neutrahno g, 
m, It is well known that recent precise determinations of the relic density Ocdm^i^ 
of non-baryonic CDM obtained by combining WMAP data with other observations of 
cosmic microwave background (CMB) anisotropics and large scale structure data, provide 
an important and often tight constraint on the CMSSM parameter space. 

It is worth remembering that the CMSSM is not the only viable and economical frame- 
work providing well defined and strongly constrained ranges of parameters. A virtue of 
the CMSSM lies in the particularly simple boundary conditions at the unification scale 
and, as a result, in a minimal number of parameters. A drawback is that the CMSSM, 
as a framework, is not broad enough to accommodate a richer structure of GUT~scale 
physics, in particular a non-minimal fiavor structure and many realistic fermion family 
patterns. One attractive and well studied model is the MSSM with boundary conditions 
at the unification scale imposed by consistency with a minimal 50(10) GUT model |^. 

The most popular approach to exploring and delimiting viable regions of the parameter 
space of the CMSSM and other SUSY models has been a usual method of evaluating the 
goodness-of-fit of points scanned using fixed grids Q. Such scans have the advantage of 
pre-determining the ranges and step size for each parameter and thus of being able to 
control exactly which points in the parameter space will be probed. On the other hand, 
the method has several strong limitations. Firstly, the number of points required scales as 

^When the CMSSM is extended to include an axionic sector, a natural candidate for the LSP and CDM 
is an axino Q], the fermionic partner of the axion. Likewise, when the CMSSM is coupled to supergravity, 
a gravitino, the fermionic partner of the graviton, arises as a possible choice for the LSP and CDM. (For 
some recent work see, e.g., |5| ^.) In contrast to the neutrahno, in both of these cases 7?-parity does not 
have to be conserved to provide a solution to the DM problem because of their very strongly suppressed 
interactions to ordinary matter. 
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, where is the number of the model's parameters and k the number of points for each 
of them. Therefore this approach becomes highly inefficient for exploring with sufficient 
resolution parameter spaces of even modest dimensionality, say > 3. Secondly, narrow 
"wedges" and similar features of parameter space can easily be missed by not setting a 
fine enough resolution (which, on the other hand, are likely to be completely unnecessary 
outside such special regions). Thirdly, extra sources of uncertainties {e.g., those due to the 
lack of precise knowledge of SM parameter values) and relevant external information {e.g., 
about the parameter range) are difficult to accommodate. 

In contrast, Bayesian statistics formalism linked to a Markov Chain Monte Carlo 
(MCMC) method of exploring a multi-dimensional parameter space offer several advan- 
tages. In the Bayesian context, the required computational power scales very favorably 
with the dimensionality of parameter space, N. The details depend on the problem un- 
der consideration, but roughly the number of points needed scales approximately linearly 
with A^. Secondly, it is straightforward to include in the analysis all sources of uncertainty 
that are present in the highly complex problem of comparing theoretical predictions with 
current collider and cosmological measurements. For instance, our imperfect knowledge of 
the relevant SM parameters has an impact on our statistical conclusions {i.e., inferences) 
about the values of the quantities of interest, in our case the CMSSM parameters. This 
source of uncertainty can fully and easily be accounted for in Bayesian statistics. Thirdly, 
probability distributions can be computed for any function of the CMSSM parameters, and 
in particular for all interesting physical observables. 

The MCMC approach has been widely used in several branches of science with much 
success and is gaining popularity in the astrophysics and cosmology community. (For some 
recent applications, see e.g. flOfl for cosmological data analysis, lTl| in the astrophysical 



context, and e.g. [12| for a general introduction to Bayesian methods.) Within the context 
of softly broken low energy SUSY, limited random (not MCMC) scans were first inter- 
pretted in the language of statistical analysis in |jl3| in the MSSM with diagonal flavor 
mass entries, and in the CMSSM in [l^, Il5|. The MCMC method was first applied to the 



CMSSM in [^] (with some modifications) and more recently in |T^, 18]. 

In the present work, we employ the MCMC algorithm to explore the parameter space 
of the CMSSM. The model is constrained with data coming from present collider data and 
cosmological observations of the CDM abundance. We compute the W boson pole mass 
Mw, sin^ Higgs and superpartner masses, BR{B — > X^'^), the anomalous magnetic 
moment of the muon {g — 2)^, BR{Bs — > and fi^/i^, and compare them with 

current experimental data. We also compute a spin-independent dark matter WIMP elastic 
scattering cross section on a free proton, <jp^, but we do not enforce upper experimental 
limits because of the uncertainties in the structure of the Galactic halo as well as in the 
values of some hadronic matrix elements entering the computation of 0"^^. We shall see 
that the current experimental limit lies just above the high-probability regions of parameter 
space. 

Our analysis goes beyond other recent works in several aspects. We include both 
experimental and theoretical uncertainties in treating relevant SM quantities. Instead of 
applying a sharp cut at experimental limits, we smear them out by incorporating theoretical 
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and experimental uncertainties. Furthermore, we improve the accuracy of CMSSM pre- 
dictions for the neutrahno rehc abundance by including previously neglected dependence 
on the remaining uncertainty in the fine structure constant. Our analysis covers larger 
values of niQ than in [^], thus allowing us to explore the focus point (FP) region |19|. We 
present inferences on the high probability regions for the CMSSM parameters and super- 
partner masses, and for the other observables listed above. We emphasize that, in some 
cases, inferences on the favored intervals of some parameters do depend on the assumed 
prior ranges, while in the other the inferred high-probability regions are fairly robust. As 
we will see, this often affects resulting conclusions about prospects for SUSY discovery at 
the LHC or in DM searches. We point out the difference between posterior probability 
(Bayesian statistics) and the quality-of-fit statistics, and emphasize that a discrepancy 
between the two methods can be only resolved with better data. 

The paper is organized as follows. In section |2| we briefly review some elements of 
Bayesian statistics and the MCMC method, and introduce our 8-dimensional parameter 
space which we explore in the following. In section |3| we describe the current collider and 
astrophysical data used in the analysis to put constraints on the CMSSM. We present 
and discuss our results in section ^. In section |5| we summarize our main findings and 
conclusions. 



2. Bayesian statistics and the CMSSM 

In this section we introduce some basic concepts of the Bayesian statistics and apply them 
to the CMSSM. 

2.1 Parameters and probabilities 

We are interested in delineating high probability regions of the CMSSM parameter space. 
We fix sign(//) = throughout (see below) and denote the remaining four free CMSSM 
parameters by the set 

6* = (mo,mi/2, Ao,tan/3). (2.1) 

In most of previous analyses, the values of the SM parameters, such as the top quark mass, 
which strongly influences some of the CMSSM predictions, have been fixed at their central 
values. However, the statistical uncertainty associated with our imperfect knowledge of the 
values of relevant SM parameters must be taken into account in order to obtain correct 
statistical conclusions on the regions of high probability for the CMSSM parameters. This 
can easily be done in a Bayesian framework by introducing a set tp of so-called "nuisance 
parameters". For the purpose of this analysis the most relevant ones are 

V' = {Mt,mb{mi>)^,aem{Mz)^,asiMz)^), (2.2) 

where Mt is the pole top quark mass, mb{mh)^^^ is the bottom quark mass at mi,, while 
aem(-^z)*^'^ and as{Mz)^^^ are the electromagnetic and the strong coupling constants at 
the Z pole mass Mz- The last three parameters are evaluated in the MS scheme. 
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The set of parameters and ip form an 8-dimensional set ij of our "basis parameters" 

V = {OA)- (2.3) 

In the foUowing, we shall specify a set ^ of several collider and cosmological observables 
which we call "derived variables", 

e=(6,6,...,em). (2.4) 

Their values depend on the CMSSM and SM parameters r] sampled in our MCMC analysis, 
^(ry). Some of the observables will be used to compare CMSSM predictions with experi- 
mental set of data d which is currently available either in the form of positive measurements 
or as limits. 

The central quantity which constitutes the basis of all probabilistic inferences is the 
posterior probability density function (pdf) p{r]\d) for the basis parameters ry. The posterior 
pdf represents our state of knowledge about the parameters r] after we have taken the data 
into consideration (hence the name). Using Bayes' theorem, the posterior pdf is given by 

pm = (2.5) 

On the rhs of Eq. (p.5|), the quantity p{d\(^), taken as a function of d for a given r/, and 
hence a given ^(r/), is called a "sampling distribution". It represents the probability of 
reproducing the data d for a fixed value of ^ (r/) . Considered instead as a function of ^ for 
fixed data d, p(d|^) is called the likelihood (where the dependence of ^ on is understood). 
The likelihood supplies the information provided by the data. In section we explain in 
detail how it is constructed in our analysis. The quantity vr(r/) denotes a prior probability 
density function (hereafter called simply a prior) which encodes our state of knowledge 
about the values of the parameters in rj before we see the data. The state of knowledge is 
then updated to the posterior via the likelihood. Finally, the quantity in the denominator is 
called evidence or model likelihood. In the context of this analysis it is only a normalization 
constant, independent of rj, and therefore will be dropped in the following. For further 



details about the terminology of Bayesian statistics, see e.g., |20, 12|. 

Since in this work we are not interested in the nuisance parameters tp themselves, at 
the end we simply marginalize over them by integrating p{r]\d) over their values. This 
procedure gives a posterior pdf for the interesting CMSSM parameters 9 which takes full 
account of the uncertainties in ip 

p{e\d) = j p{e,i)\d)d^i). (2.6) 

Note that all pdf's should normally be normalized so that the total probability is unity. 
However, for the parameter estimation procedure presented here only the relative posterior 
pdf's are relevant, and in the following we shall plot pdf's normalized in such a way that 
their maximum value is one. In practice, we will present pdf's for only one or two variables 
at the time, with the remaining ones integrated over. We will introduce and discuss several 
of them in Section Q where we present our results. 
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The purpose of the MCMC exploration of parameter space is to obtain a series of points 
(called a "chain"), whose density distribution is proportional to the posterior pdf on the rhs 
of Eq. ( |2.5| ) . Further details about the MCMC procedure are given in appendix From 
the samples in the chain, it is straightforward to obtain all pdf's of interest by plotting 
histograms of the number of samples as a function of the parameter values that one wants 
to examine. In particular, one does not need to carry out the marginalization integral in 
Eq. ( p. 61 ) explicitly. It is sufficient to ignore the coordinates of the samples in parameter 
space along the marginalized directions. This is one more major advantage of the MCMC 
method. 

Another useful feature is that from the posterior pdi p{ri\d) we can obtain the posterior 
pdf for any function / of basis parameters, using the fact that 

p{f, V\d) = p{f\v, d)p{v\d) = difiv) - f)p{v\d), (2.7) 

where 5{x) denotes the delta-function. Therefore for every sample in the MC chain, one 
simply computes /(r/) and the resulting density of points in the (/, tj) space is proportional 
to the posterior pdi p{f , r]\d) . From this pdf one can then obtain by marginalization the 
pdf for any subset of {f-,Tj). In particular, if we take /(??) = ^(j?), we can then investigate 
probability distributions for any combination of the basis and the derived parameters, as 
well as their correlations. This is investigated in section ^. 

Before we can proceed to delineating high probability regions in the CMSSM parameter 
space, first we must choose our priors and specify the likelihood. This is the subject of the 
next two sections. 

2.2 The choice of prior probabilities 



It is clear from the rhs of Eq. (|2.5D that in the Bayesian approach the first step is to 
specify the functional form of the prior pdf. This is equivalent to assigning a probability 
measure to parameter space. The principle of indifference states that one should assign 
equal probabilities to equal states of knowledge before seeing the data. In our case, the 
basis parameters are location parameters over which it is appropriate to set a flat prior 



■K{ri) 



const for r/min < l] < ?7max 

otherwise. 



where the constant is determined by the requirement that the prior integrates to probability 
one. Since we assume no correlation of priors between the SM and the CMSSM parameters, 
the joint prior can be written as 

7r(r/) = 7r(6l)7r(V'). (2.9) 

Flat priors are thus characterized by their ranges [^/mim ??max] • One alternative possibility 
would be to employ a "naturalness" prior which gives more weight to points exhibiting less 
fine-tuning Q. 

Our SUSY and SM parameter priors are summarized in Table ||. In this work we 
consider two initial ranges of CMSSM parameters. In one, which we call a "2 TeV range". 
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CMSSM parameters 9 



"2 TeV range" 



50 GeV < mo < 2 TeV 
50 GeV < mi/2 < 2 TeV 
Uol < 5 TeV 



"4 TeV range" 



50 GeV < mo < 4 TeV 
50 GeV < mi/2 < 4 TeV 
\Ao\ < 7 TeV 
2 < tan /3 < 62 



SM (nuisance) parameters ijj 



160 GeV <Mt< 190 GeV 
4 GeV < m5(mb)^ < 5 GeV 
127.5 < l/a^^iMzY^ < 128.5 

0.10 < OsiMzY^ < 0.13 



Table 1: Initial ranges for our basis parameters 77 ~ {9, 1}}), with flat prior probability distributions 
assumed. 



we assume 50 GeV < mo,mi/2 < 2 TeV and l^o] < 5 TeV. This choice is motivated by 
an expected LHC reach in exploring superpartner mass and by a general "naturalness" 
argument of SUSY mass parameters to preferably lie within 0(1 TeV). In the other case, 
called a "4 TeV range", we assume 50 GeV < mo,mi/2 < 4 TeV and |^o| < 7 TeV, 
which goes far beyond the LHC reach. (The larger range will include the focus point 
region, along with various uncertainties involved. We will discuss this point later.) We will 
compare our findings for both ranges in order to see to what extent statistical conclusions 
depend on our preconceived expectation that a SUSY signal might lie within reach of the 
LHC (represented by the "2 TeV range"). Such a sensitivity test is essential to establish 
the extent to which inferences depend on the initial range one chooses, i.e. on the prior. 
The lower bounds on mo and mi 12 come from the negative results of sparticle searches. 
We allow a rather generous range for Aq, in part to see to what extent this choice would 
allow one to reduce [21| the impact of the cosmological constraint, which in the usually 
explored case of Ao = is very tight. For both sets we further assume 2 < tan (3 < 62. 
The lower bound comes from negative Higgs searches p^ . Very large values of tan /? > 60 
are in conflict with theoretical considerations, e.g. they would make it extremely diflicult 
to achieve radiative electroweak symmetry breaking [23|. Furthermore, at such large tan/? 
large uncertainties arise in the computation of the SUSY spectrum, leading to unreliable 
predictions. On the other hand, since the SM nuisance parameters are well measured, it 
turns out that their prior ranges are irrelevant for the outcome of the analysis. 

Before closing this section, we comment that the necessity of choosing priors is often 
(incorrectly) regarded as a limitation to the "objectivity" of the Bayesian approach. This 
can be easily dispelled by noting that two scientists in the same state of knowledge before 
seeing the data (i.e., who have assumed the same priors) will necessarily reach the same 
conclusions. When the choice of priors makes a difference in drawing the final inference 
(given by the posterior pdf), this is a "health warning" that the data is not informative 
enough, e.g. the likelihood is not sufficiently peaked to override the assumed prior distri- 
bution. In this case, the inference on the parameters must rely either on external relevant 
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information in the form of the prior {e.g., a theoretical "naturalness criterion" investigated 
in Ref. or of better and more constraining data. In the present study, this will be 

the case for the parameter mo, as we discuss in section 4.1. 



3. Collider and cosmological observables 

In this section, we first define the likelihood function for the CMSSM. Next, we introduce 
the data that we use for the nuisance parameters, and then proceed to describe electroweak 
and dark matter observables. We give details about their calculation, the theoretical un- 
certainties involved and the experimental errors in their determination. 

3.1 Constructing the likelihood for the CMSSM 

The likelihood is a key element of our analysis. It encodes the information from the ob- 
servational data and therefore particular care must be taken in constructing it. In the 
Bayesian framework it would be easy to incorporate the full likelihood functions from 
various experimental measurements if they were available. However, even though the ac- 
tual measurements contain much more useful information, most measurements in particle 
physics experiments are presented only by the mean and the standard deviation, while 
upper or lower exclusion bounds are usually given in terms of the 95% exclusion CL. 

Uncertainties in the observable quantities can be split into two categories. The first is 
an experimental uncertainty, the second is a theoretical one and is a consequence of making 
some approximations {e.g., neglecting higher order loop corrections), a limited numerical 
precision in the code, etc. In practice, the theoretical uncertainty can be modelled in a 
Bayesian context by considering that our mapping from the basis parameters ry to the 
derived quantities ^ is imperfect: instead of an "exact" mapping ^(ry) we actually have 
only an imperfect version ^(77) which suffers from the sort of uncertainties outlined above. 
The likelihood p{d\(,) introduced in Eq. ( |2.5| ) can then be written as 

p{d\0 = I p{d\Op{md"'i (3.1) 

where we have integrated out the true (and unknown) mapping. The pdf p{C\0 encodes 
the estimated uncertainty of our mapping. Usually we only have (at best) an estimate 
of the theoretical errors, which means that we only have information on the scale of the 
associated uncertainty. This is described by a multi-normal distribution of a general form 

where C is an m x m covariance matrix describing the error of the mapping (m being 
the number of elements in ^) and \C\ denotes its determinant. If one assumes that the 
theoretical errors Tj {i = l,...,m) for the different quantities are uncorrelated then C 
is diagonal, C = diag (ti, . . . , r^). Furthermore, if the likelihood p{d\S,) is also a multi- 
dimensional Gaussian function with diagonal covariance matrix, then we have 

Pidli) = exp (-^(d - i)D-\d - , (3.3) 
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where D = diag (cti, . . . , 0"^) and cij denotes the experimental standard error. In this 
simple case Eq. ( |3.lD reduces to the usual rule of adding theoretical and experimental 
errors in quadrature for each derived observable, i.e., the total error in each direction is 
~ Y^cf + . This familiar result is a special case of the more general treatment given 
above, which shows that in the Bayesian framework all sources of uncertainty specified in 
the model can be fully taken into account in a systematic way. 

In this work we model the likelihood of all observables for which there exists a positive 
measurement by an uncorrelated, multi-dimensional Gaussian function. The experimental 
means and standard deviations for SM quantities are summarized in Table and will be 



discussed in more detail in section |3.2| . Since the SM quantities are at the same time input 
quantities of our analysis, there are no theoretical uncertainties associated with them. In 
Table |3| we display the experimental and theoretical errors for derived cosmological and 



collider quantities, see section for a further discussion. 

For the quantities for which only lower or upper bounds are available {e.g., superpartner 
masses or BR{Bs — > fi^fi~)), a usual procedure is to simply discard points in the parameter 
space for which such limits are violated at some confidence level {e.g., la or 95%), essentially 
using a step function as a likelihood. This procedure is not totally rigorous since it does not 
take into account the amount of uncertainty associated with the theoretical error (denoted 
by r). In our analysis instead we use Eq. (|3.1| ) to incorporate our estimate of the theoretical 
uncertainty of the mapping in the limits that we use, as explained below. 

As an illustration, let us consider a 1-dimensional case involving only one observable 
^. A lower bound on the (exact) mapping ^ can be described by the following likelihood 
function (replacing d — > (d, ^Hm)): 



P(0-,6im|6 



1 



27ro- 
1 



2iTa 



exp 



(e-ei 



limy 



2^2 



for ^ > Ciim, 
for i < ^lim. 



(3.4) 



where a Gaussian function has been used to model the drop of the likelihood function 
below the experimental bound ^um- As explained above, the theoretical uncertainty of the 
mapping is described by a 1-dimensional Gaussian of standard deviation r (see Eq. (| 
for m = 1). Then from the integral in Eq. ( |3.1| ) we obtain 



where ^ 



exp 



S,{9) and we have defined 
T Vo-2 + r2 



(6im-e)^ 

2(a2 + T2) 



Z{tlim) = 



[1 - Z(tHm)] + Z 



dx exp{—x'^ /2). 



(3.5) 



(3.6) 



*lin 



This form of the likelihood encodes the uncertainty associated with our imperfect 
mapping, as described by r. The effect of including the theoretical uncertainty is to smear 
out the drop of the likelihood function: the scale of the drop goes from a to Vo"^ + 
Unfortunately, experimental bounds are usually given in the form of the lower or upper 95% 
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Figure 1: An illustration of the likelihood function p(cr, ^uml^', t) used for a quantity for which only 
a lower bound is available, including the theoretical uncertainty t and setting the experimental error 
a = 0. The dashed green line is the sharp 95% CL bound {£,\im = l), the solid red curve includes 
the theoretical uncertainty of r = 0.05 which smears out the limit. 

confidence limit only, without the possibility of deriving a form of the likelihood analogous 
to ( p. 51 ). In the absence of fuller information about the experimental value of a, we can take 
in ( |3.5| ) the limit a <^ t and simply use the 95% CL bound as ^um- This procedure leads 
to the likelihood function plotted in Fig. 0, showing that the inclusion of the theoretical 
uncertainty smears the 95% bound over the scale r. This procedure is more conservative 
than the usual method of simply rejecting all points below ^nm. 

We use ( |3.5| ) as the likelihood function for all of the upper/lower bounds listed in 
Table |^ along with our estimates of the corresponding theoretical uncertainty. We also 
discard, i.e., assign zero likelihood to all unphysical points in parameter space, i.e., those 
for which any of the masses becomes tachyonic or the conditions of EWSB are not satisfied. 
We do the same for the cases where the lightest neutralino is not the LSP. 

3.2 Inputs for SM quantities 

As explained in section |^ the uncertainties coming from experimental errors in the SM pa- 
rameters are incorporated through four nuisance parameters: Mt, mh{mi))^'^^ , aemiMz)^^ 
and as{Mz)^^^ . These quantities and their uncertainties are summarized in Table |2|. 

Note that for the running bottom quark mass niiy{rrn,)^'^ ^ the range adopted in Table § 
is rather conservative. It arises from combining measurements from bb systems, b-flavored 
hadrons and high-energy processes. More recently, a much more precise determination of 
mi,{mb)^^'^ = 4.19 lb 0.06 GeV has been obtained from a renormalization group improved 
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IVTpa n v?i 1 n p T In rprta i n tv 
/i cr (exper.) 


Ref. 


Mt _ 


172.7 GeV 2.9 GeV 
4.24 GeV 0.11 GeV 
0.1186 0.002 
127.958 0.048 


1 
1 
1 
1 


21 

25| 
26| 
27| 



Table 2: Experimental mean and standard deviation a for the nuisance parameters used in the 
analysis. For all these quantities we use a Gaussian likelihood function with the mean /i and the 
standard deviation a. 



Derived observable 


Mean value Uncertainties 

/i cr (exper.) r (theor.) 


Ref. 


Mw 
siv? 

BR(i^ ^ Xs-i) X 10"^ 


80.425 GeV 34 MeV 13 MeV 
0.23150 16 X 10-^ 25 x 10"^ 
25.2 9.2 1 
3.39 0.30 0.30 
0.119 0.009 Q.lVL^h'^ 


|27 

w^ 

[30 
|32 
Th 


31] 
IS work 



Table 3: Summary of derived observables used in the analysis for which positive measurements 
have been made. As explained in the text, for each quantity we use a likelihood function with 
mean jj, and standard deviation s = Vcr^ + t^, where a is the experimental uncertainty and r our 
estimate of the theoretical uncertainty. 



sum rule analysis p8|. As for the fine structure constant, the dominant error comes from 



the hadronic contribution Aahad = 0.02758 it 0.00035 |27|. Note that, because of this 
uncertainty OemiMz)'^^^ is not as precisely known at the scale Mz as at zero momentum 
transfer. Not included in the list of nuisance parameters are the pole mass of the Z boson 
Mz = 91.1876(21) GeV, as weh as the Fermi constant Gf = 1.16637(1) x 10"^ GeV'^ 
As one can see, they are now known with very high precision and we fix them at their central 
values. Including them as nuisance parameters would not have any appreciable effect. 

3.3 Derived observables 



We now present the derived observables which in section 2.1 were denoted by a general 
symbol ^. In our analysis they are computed in terms of the CMSSM parameters: mg, 
"^1/21 ^0) tan/3, as well as the nuisance parameters discussed above. 

In Table ^ we summarize the derived variables for which positive measurements have 
been made while in Table § we list the ones for which currently only experimental bounds 
exist. Lower bounds on sfermion masses have been obtained in the context of the MSSM 
at LEP-II and Tevatron Run II. Below we comment on some of the entries and on our 
procedure for their computation. Generally, to calculate Higgs and SUSY mass spectra 



we use the package SOFTSUSY vl.9 [g9[ which employs 2-loop RGEs for couphngs (both 



gauge and Yukawa) as well as for gaugino and sfermion masses. 

W gauge boson mass Intrinsic theoretical uncertainties coming from higher loop effects 
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Derived observable 







TTT, 






T TT 


rrih 




T T 
IjLi 


Ch 




TTT 






T T 






T T 

LiLi 


"'xf 






mg^ 




T T 






LL 






LL 






LL 






LL 






LL 


rriq 




LL 


rUg 




LL 



Constraints 



Clim 



r (theor.) 



WIMP mass dependent 
1.5 X 10""^ 

114.4 GeV (91.0 GeV) 
50 GeV 

103.5 GeV (92.4 GeV) 
100 GeV (73 GeV) 

95 GeV (73 GeV) 
87 GeV (73 GeV) 

94 GeV (43 GeV) 

95 GeV (65 GeV) 
95 GeV (59 GeV) 
318 GeV 

233 GeV 



- 100% 
14% 
3 GeV 



Ref. 



135 

Wi 

i 

l34 



m 



135 

i 
1 

l3£ 



|35 



1 35 



m) 



Table 4: Summary of derived observables for which only limits exist, with UL — upper limit, 
LL = lower limit (at ^um (3.4), 95% CL, unless otherwise stated). The experimental limit on the 



spin-independent cross section for a WIMP elastic scattering on a free proton a^^ is not included 
in the likelihood, as explained in the text. Since the precise form of the likelihood is not available, 
we use the conservative procedure of including at least an estimated theoretical uncertainty r. The 



likelihood is then given by Eq. (3.5), in the limit where ct <C t. The value in parenthesis indicates 
the more conservative bound, see the text for details. Note that theoretical errors for scalar masses 
are probably much larger in the focus point region, as discussed in the text. 



to the W boson pole mass are estimated to be t{Mw) = 13 MeV [40, 41 1. The parametric 
uncertainties are dominated by the experimental error of the top quark mass and by the 
hadronic contribution to the shift of the fine structure constant. Both uncertainties are 
fully taken into account (via Aq/j^^) by marginalizing over Mt as a nuisance parameter. 
We compute Mw in the DR scheme (from gauge couplings relations), including full 1-loop 
contributions [40|. 

Effective leptonic weak mixing angle The effective leptonic weak mixing angle sin^ 9^^ 
receives SM and SUSY contributions. In our computation we include a full 1-loop SM and 
full 1-loop universal Z-vertex supersymmetric corrections EH] . The net contribution of the 



non-universal corrections is negligible |42]. Intrinsic theoretical uncertainties come from 
higher-loop effects, which induce an uncertainty^ taken to be t{sw? = 25 x 10^^. 

Note that we treat aem(-Wz)*^"^ as a nuisance parameter but take Mw and sin^ 6^^ as 
derived quantities. This is because we use a parametrization in which the former (along 
with Mz and Gp) are taken as inputs from which one can compute the latter and compare 
with experiment. 



^In Ref. a smaller uncertainty is used (i.e., r(sin'^ 9^g) = 12 x 10 ^) as a result of including leading 
two-loop supersymmetric corrections. 
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The anomalous magnetic moment of the muon The final measurement of the 
anomalous magnetic moment of the muon, = {g — 2)^, by the Brookhaven E821 ex- 
periment a^^^ = (11659208.0 it 5.8) x lO"^'^ remains in an apparent disagreement with 
SM predictions. The current SM theoretical value, based on e+e~ low energy data, is |31| 



af^ = 11659182.8 ± 6.3had ± 3.5lbl ± 0.3qed+ew) x 10"^°. The discrepancy between the 
two, if confirmed, could be attributed to an additional contribution from loops involving 
superpartners, 

JflSUSY ^ ^exp _ ^SM ^ (25.2 ± 9.2) X 10-1°, (3.7) 

where the independent errors have been added in quadrature. 

In calculating ^a^^^^ we have taken into account the full one-loop SM+SUSY con- 



tributions [43 1 and several two-loop corrections. The first class of two-loop corrections 



comprises the leading one-loop diagrams with a photon in the second loop |44|. The sec 



ond class comprises diagrams with a closed loop of SM fermions or scalar fermions 
The last class comes from diagrams containing a closed chargino-neutralino loop |^^. As 
a consequence, in the CMSSM parameter space the intrinsic uncertainties are estimated 
to be T((5a^^^'^) = 1 x 10"^'^ [47]. The sign of the SUSY contribution to is the same 
as the sign of fi. Since the former is positive, throughout this analysis we have assumed 
sign(;u) = +1. 

Finally, we note that, the SM value evaluated using r-data is = (11659201.8 + 6.3, 
leading to only a 0.7o" deviation [^]. In light of this we do not feel it is justified to apply 
a SM prediction based on combining the two sets of data. Secondly, there appears to be 
potentially some discrepancy among different sets of e^e~ which again may put the claimed 
deviation ( ^ ) from the SM into question. Because of these outstanding problems, in 
subsection we will perform a separate analysis where we exclude from the likelihood. 

BR{B — y Xg'-f) The current experimental world average value by the Heavy Flavour 
Averaging Group (HFAG) is [|2[ 

BR{B ^ X,7)exp = (3.391°:^°) x 10-^ (3.8) 

which agrees rather well with the full NLO prediction of the SM [^] 

BR{B X,7)sM = (3.70 ± 0.30) x 10"^. (3.9) 

This clearly imposes an important constraint on any additional contributions. In the 
CMSSM, the assumed universality of soft mass terms leads to a particularly simple flavor 
structure in which no additional sources of flavor mixings beyond those due to the CKM 
matrix are present, the framework known as minimal flavor violation (MFV). In this case 
SUSY contributions arise from a loop involving the top quark and the charged Higgs boson 
and one of the stop-chargino exchange. In a more general flavor mixing scenario additional 
one-loop contributions arise due to gluino (or neutralino) and down-type squark exchange. 

We compute SUSY contribution to BR{B ^sT) following the procedure outlined 



in Refs. [5C, ^ where, in addition to full leading log corrections, large tan /3-enhanced 



terms arising from corrections coming from beyond the leading order (BLO) have been 
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included.'^ Furthermore, combined experimental constraints from BR{B — > Xg")) and 
BR{B — > Xsl'^l~) imply that, baring highly non-standard scenarios |5^, the sign of the 
total amplitude for the decay B —>■ Xg'y has to be the same as that of the SM 

We assume no significant additional theoretical uncertainty beyond that coming from 
the SM calculation. This can in part be justified by the fact that the overall SUSY 
contribution to BR{B — > Xg^) has to be small, as noted above. More importantly, we 
estimate that most of the uncertainties due to SUSY contributions are due to remain- 
ing uncertainties in the values of SM parameters, especially the top and bottom masses 
which are accounted for in their treatment as nuisance parameters. Therefore we take 
t{BR{B Xg-i)) = 0.30 X 10-^. 

Cosmological constraints on the dark matter density A combination of the recent 



WMAP data |56] on the CMB anisotropics with other cosmological observations, such as 
measurements of the matter power spectrum, leads to tight constraints on the cold dark 
matter (CDM) relic abundance. The exact numerical value depends on a number of as- 
sumptions about the underlying cosmology [e.g., the geometry of space, the adiabaticity of 
initial conditions and a power-law, feature-free primordial power spectrum). It is impor- 
tant to keep in mind that relaxing some or all of these assumptions can considerably weaken 
the constraints. Furthermore, different combinations of data sets also lead to somewhat dif- 
ferent values for the relic abundance. For definiteness, we have performed a re-analysis of 



CMB data from the first year of WMAP observations [56|, as well as from CBI |57], VSA |58| 
and ACBAR [^] , and combined them with the real-space power spectrum of galaxies from 
the SLOAN galaxy redshift survey (SDSS) |Q, restricted to scales over which the fluctu- 
ations are assumed to be in the linear regime, i.e. for k < 0.1/i"^Mpc, with the Hubble 
Space Telescope measurement of the Hubble parameter Hq {Hq = 100 /ikm/sec/Mpc), 
and with the latest supernovae observations data [^]. Assuming a flat Universe ACDM 
cosmology, the resulting constraint on the dark matter relic abundance - after marginaliz- 
ing over all other relevant cosmological parameters - is well approximated by a Gaussian 
function with a mean and standard deviation given by"^ 

J^cdm/i^ = 0.119 ±0.009. (3.10) 

We make use of this constraint and we assume that all the CDM is made up of stable 
neutralinos, but we also enlarge the error to include a theoretical uncertainty in the com- 
putation of r^cDM^^ (see below). 

A precise determination of the neutralino relic abundance fi-^/i^ requires an accurate 
treatment of the neutralino pair annihilation and coannihilation cross sections into all 
SM particle final states. We employ exact expressions for neutralino pair annihilation 
processes into all allowed final-state channels which have been computed in [^] and which 



^An analogous and updated BLO-level analysis in the case of general flavor mixing has been performed 
in Refs. [ p2| , p^ , p4[ and shows in general a much larger difference between LO and BLO results than in 
MFV. 

1, 4-„„l ,,„l..„ ,„ ,,l.„u-n,, J.fr„„ 4- c , o u2 _ n no 



"The central value is slightly different from QcDMh'^ = 0.113io:oo9 obtained by the WMAP team in §6 



because of the different combination of data employed. 
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are valid both near and further away from resonances and thresholds. We further include 
the neutralino coannihilation with the lightest chargino and next -to-lightest neutralino 1 64 1 
and with the lighter stau [|65| ] with similar precision. We include all coannihilation channels, 
including coannihilation with light stops [36|. We compute Q^h'^ by solving the Boltzmann 
equation numerically as in [S7|. 

As is well known from fixed-grid scans, in the CMSSM the values of J^^^^ typically 
exceeds the range given in (3.10), except in some rather special regions of the parameter 
space. Firstly, at fairly small mi/2 and mo, the neutralino annihilation channel into SM 
fermion/boson pairs via t/u-channel exchange of a superpartner opens up only a narrow 
band consistent with J^cdm^^ (the "bulk" region) which however is largely excluded by 
a lower bound on the Higgs mass and/or by a strong upper bound on allowed SUSY 
contribution to BR{B — > Xg^). Secondly, if the LSP and the next-to-LSP (NLSP) are 
closely degenerate in mass, then the LSP coannihilations with NLSPs near freeze-out may 
reduce the LSP relic density considerably. In the CMSSM, efficient coannihilation takes 
place along the boundary dividing neutralino and lightest stau ri LSP regions. In addition, 
at large Aq there are limited cases where the lightest stop is the NLSP and is almost 
degenerate in mass with the LSP |6^]. Thirdly, in some cases neutralino annihilation can 
be enhanced via the process XX ~^ ff involving an s-channel Higgs or Z boson exchange. 
At large tan /? ~ 50 the most significant effect comes from the CP-odd Higgs A resonance 
(around niA — 2m;^,, where is the lightest neutralino mass), since the A couplings 
to down-type fermions are enhanced at large tan f3, and because, in contrast to heavy 
scalar Higgs exchange, the process is not p-wave suppressed. Finally, at very large rriQ of 
a few TeV the rapidly decreasing fj? (from large positive to negative values) causes the 
higgsino component of the LSP to increase which, in a narrow focus point (FP) region (of 
still positive /i^) allows Qy^h'^ to pass through the favored range (|3lO| ), before becommg 
too small. 

The theoretical uncertainty involved in computing Qy^h'^ varies greatly depending on 
the case. Errors come from computing the Higgs and superpartner mass spectrum at 
finite (two-loop) order in RGEs, a scale dependence, finite numerical accuracy in solving 
the Boltzmann equation, some residual uncertainties in computing the gauge couplings, 
as well as potentially much larger errors in computing top and bottom Yukawa couplings 
and Higgs widths. The choice of the scale at which one minimizes the effective potential 
has a minor effect [17]. A numerical algorithm of solving the Boltzmann equation is very 
accurate and in the bulk region has an estimated error of a few per cent ||68| , |63| , which is 
comparable with the observational error in ( 3.10D . 

In the bulk region where fi^/i^ primarily depends on the t/ii-channel exchange of slep- 
tons, uncertainties in the calculation of the SUSY spectrum are 0{1%) (for moderate values 
of the CMSSM mass parameters, even at large tan/3) and therefore under control [^]. On 
the other hand, the accuracy is much poorer in the regions of special cosmological in- 
terest that have been mentioned above. In the coannihilation regions Jl^^/i^ is sensitive 
to e"^™'/'^, where Am = rriNLSP ~ difference between the mass of the NLSP 

"iNLSP and m^. Even a small variation of 0(1%) in Am can lead to ~ 30% variations in 
[0]' This sort of error is inherent in current determinations of the mass of ti which 
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is the NLSP in large parts of the {1^1/2, nT-o) plane. Larger uncertainties arise if the NLSP 
is the ti since the shift in its mass 10% can lead to order-of-magnitude discrepancies in 
the prediction of the relic density . In the CP-odd Higgs resonance region one finds a 
strong suppression of ^x^'^ broad ranges of mi/2 and niQ. In the CMSSM this occurs 
for large tan /3 where effects of the bottom Yukawa coupling hh on the RGE running are 
important and reduce niA compared to the low tan/? regime. The relic density is further 
suppressed by the enhancement of the coupling Abb which is proportional to hb- Special 
attention must also be paid to the computation of the Higgs boson width which receives 
sizable radiative corrections. We include one-loop QCD corrections ||7^ as well as those 
due to (QCD corrected) Yukawa vertices. Still, unknown two-loop corrections to hb may 



cause an uncertainty of up to 30% [71|. In the FP region, a determination of /i is strongly 
sensitive to two-loop corrections proportional to ht, whose computation requires special 
care, and can lead to ~ 100% errors in the computation of J^^^^ lOl- 

In addition to the theoretical errors discussed above, the value of J^x^^ depends on 
the top and bottom masses, as well as aem(Afz)^'^- In particular, we have found that 
the seemingly small experimental error in aem{Mz)'^^^ leads, indirectly, via its effect on 
the SM gauge couplings gi^2i to a variation in il^/i^ of order 10%. All of these effects are 
accounted for by our use of the nuisance parameters. 

The error in the quantity ^x^"^ thus be very sensitive to which (co) annihilation 
process is most efficient. This makes it difficult to evaluate the theoretical uncertainty in 
O^/i^. Given the above discussion we estimate r(r2^/i^) = 10% in the bulk of the parameter 
space although we are aware that in the FP region the error is almost certainly much larger. 
We add this error in quadrature to the observational error cr{QcDMh'^) = 0.009 in ( 3.1C| ) 
and obtain the Gaussian with the following mean and standard deviation. 



^xh^ = 0.119 ± J(0.009)2 + (0.1 17^/12)2 (3.11) 



= 0.119 ± 0.009 yi + 1.75 (S7^/i2/0.119)2. (3.12) 

Note that the theoretical uncertainty is of the same order as the uncertainty from current 
cosmological determinations of JIcdm^^- 



We now comment on some of the derived variables for which the 95% CL experimental 
upper/lower limits have been presented in Table We start, however, with the only 
observable which is not included in the likelihood. 

Direct detection of dark matter Assuming the Galactic DM halo is mostly made up of 
neutralinos, it may be possible to directly detect them via their elastic scatterings off nuclei, 
or indirectly via their annihilation products. Direct DM detection in SUSY frameworks 



has been investigated by many authors |72]. The CDMS Collaboration (CDMS-II) has 



recently improved their previous upper limit on the spin-independent dark matter WIMP 
elastic scattering cross section on a free proton, cXp^, down to some 2 x 10~^pb (at low 
WIMP mass) |7^. Comparable upper limits have also been set up by the Edelweiss-I ||7^ 



and ZEPLIN-I |75| experiments 
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However, several questions regarding the properties of the DM halo {e.g., the existence 
of clumps of DM and the value of the local halo mass density) remain unsettled. Recent 
numerical N-body simulations of large structure formation have revealed a large number of 
overdense regions surviving until today |7^. It is possible that an improved sensitivity of 
the simulations will reveal further, and finer, dumpiness. The clumps may contain a sizable 
fraction of the total dark matter halo, of the order of 10%. It is therefore not unlikely that 
locally (at the Earth's location) the DM density may be significantly different from the 
usually assumed average value of = 0.3 GeV/cm^. 

It is worth remembering that, in translating null experimental results for an elastic 
WIMP-target cross section into upper limits on <Tp^ , not only the local DM density enters 
but also a WIMP velocity distribution. This is usually taken to be Maxwellian with a peak 
at 220 km/sec with an estimated error of between 20 and 50 km/sec. While this leads 
to an additional uncertainty in 0"^^, it is actually tiny compared with the uncertainty of 
the actual WIMP density at the position of the Earth caused by cuspiness. In view of the 
above discussion, in our opinion upper experimental limits on 0"^^ should not be used to 
constrain supersymmetric parameters with the same degree of reliability as from collider 
searches. 

In computing 0"^^ in the CMSSM, we include full supersymmetric contributions which 



af = (3.13) 



have been derived by several groups [77, 75, 8C, 81|. 0"^^ can be expressed as 

4 

TT 

where fip = mpm^/{mp + m^) is the reduced mass of the WIMP-proton system. (For 
spin-independent interactions of neutralinos, and more generally Majorana WIMPs, there 
is no need to consider an analogous quantity on a free neutron since crf^ = (Tp^ .) 



The coefficients fp^n can be expressed as |£7| 



fp firq , 2 (p) fq 

u,d,s ^ q=c,b,t ^ 



where f!^^ = 1 — '^q=uds fxq^ ^^'^ nuclear form factors f!j?^ are defined via {plmgqqlp) = 
rripflf^ {q = u, d, s), and analogously for the neutron. The masses and ratios Bg = {p\qq\p) 
of light constituent quarks in a nucleon come with some uncertainties. For definiteness, 
we adopt the set of input parameter given in [^] and assume ijiu/md = 0.553 ± 0.043, 
ms/nid = 18.9 ± 0.8, and Bd/B^ = 0.73 ± 0.02, as well as 

/^P) = 0.020 ± 0.004, = 0.026 ± 0.005, f^] = 0.118 ± 0.062, (3.14) 

f^^ = 0.014 ± 0.003, f!p] = 0.036 ± 0.008, f^J = 0.118 ± 0.062, (3.15) 

and for the parton density functions we employ the CTEQ6L set Is^ evaluated at the 



QCD scale defined by the averaged squark mass and neutralino mass Q y ~ "^x- 

The nuclear form factors f^g^^ come with a large error and these are the ones that provide 
the dominant contribution to 0"^^. This has a rather significant impact on the size of 
(Tp^ [83|, unless tan/3 is very small. 
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BR{Bs — >■ ) The latest 95% CL experimental upper limits from the D0 Run II 

and CDF Run II experiments at Fermilab are, respectively, 



BR{Bs /u+^") < 2.0 X 10"^ (CDF) ||, (3.16) 
BR{Bs /u+^i") < 3.7 X 10"^ (D0) M]. (3.17) 



A combined 95% CL limit is BR{Bs fi^i-i ) < 1.5 x 10 Ultimately, assuming the 

integrated luminosity of 8fb^^, a combined CDF and D0 limit is expected to reach some 



2 X 10 ^ [88| which is significantly above the SM prediction 

BR{Bs ^+Ai")sM = (3.42 ± 0.54) x 10"^ (3.18) 

We compute the SUSY contribution to BR{Bs — > ^~^//~) by following a full one- 
loop calculation of Ref. [^] which assumes MFV. Furthermore we include tan j3 enhanced 
corrections to the bottom quark mass |^l[. The parametric uncertainties are associated 
with an error in the decay constant /b^,, which arises from lattice calculations, and an error 



in the bottom quark mass |92]. The latter is accounted for by the MC procedure, while the 
former is of order 10%. Unknown higher order corrections are of order 10% ||9^ . Therefore 
we use a total theoretical error of 14%, obtained by adding the above uncertainties in 
quadrature. 

The lightest MSSM Higgs boson mass A final LEP-II lower bound on the SM Higgs 
mass is ?n.//g^ > 114.4 GeV (95 % CL) |2^. The bound applies to the lightest Higgs boson 
h in the MSSM if its coupling to the Z boson is SM-like, i.e. if = Ozzh/dzzHsm ~ 
sin^(/3 — a) ~ 1. This occurs in the decoupling regime where ^ mz- For arbitrary 
values of m^, the LEP-II Collaboration has set 95% CL bounds on mh and rrij^ as a 
function of [22|, with the lower bound of m/j > 91 GeV for ~ tua and <^ 1 [22|. 



In this low-mass region we use a cubic spline to interpolate between some selected points in 
nifi and derive the corresponding 95% CL bound, which is then smeared with a theoretical 
uncertainty r of 3%. The intrinsic theoretical error in computing ruh, after taking into 
account effects of the renormalization scale dependence, in the CMSSM has been estimated 



to be r(m/i) = 3 GeV [^, 95]. The parametric uncertainty coming from the errors in top 
quark mass and the strong coupling constant are accounted for by taking them as nuisance 
parameters. 

As mentioned above, we compute the lightest Higgs mass following the SOFTSUSY vl.9 
package , where full one-loop and leading two-loop corrections and two-loop effects on 



the EWSB conditions are included. 

Superpartner masses Below we comment on some limits on superpartner masses. Since 
currently there is no information available on the likelihood function for sparticle masses, 
we make use of the experimental error of 95% CL and of the likelihood given in Eq. (p.5|). 

The parametric uncertainties in the sparticle masses coming from the SM variables are 
accounted for via the nuisance parameters. The authors of Ref. have argued that the 
theoretical uncertainties in the computation of SUSY masses are of order 0(1%) except 
in some special regions of the parameter space (such as the FP region) which require a 
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separate treatment. Conservatively we take r = 5% for each computed superpartner mass 
and use it to smear out the 95% CL experimental lower limit on its mass, as explained in 



section 3.1 



Neutralino LSP mass In the context of the CMSSM, LEP-II provides an absolute 
lower bound on the mass of the lightest neutralino LSP xi (in this work denoted simply 

by x) H 

(3.19) 



> 50 GeV. 



Chargino mass Chargino mass has been excluded up to m ± > 103.5 GeV [34|, provided 
that mj> > 300 GeV, where mj> stands for the lightest sneutrino mass which in the CMSSM 
is Dr- However, when the mass difference ?7i^± — < 3 GeV, as in the FP region, then 



the bound is relaxed to m^± > 92.4 GeV [34]. The latter bound is also applied when 
m,y < 300 GeV. 



Slepton masses Combined slepton mass limits have been obtained by the LEP SUSY 
working group [^]. The overall limits are 



me- > 100 GeV, 



m^jj > 95 GeV, 



ruf^ > 87 GeV. 



(3.20) 



These limits are valid provided that mr — rriy > 10 GeV. Otherwise we apply the more 



conservative bound m^-^ > 73 GeV 36]. 

For the sneutrino, in the context of the CMSSM, the DELPHI collaboration has ob- 
tained the following limit J^] 

my > 94 GeV, (3.21) 



provided mc — > 10 GeV. Otherwise we apply m/j > 43 GeV ]p6 ] 



Squark masses Combined limits on the lightest stop and sbottom masses from LEP-II 
are JQ] 

m^^ -^^ > 95 GeV, (3.22) 

provided m^^ ^ GeV. Otherwise we apply the more conservative limits obtained 

by the ALEPH Collaboration ]p5| 



(3.23) 



> 65 GeV and m^^ > 59 GeV. 



Finally, for the two first generations, the D0 Run II Collaboration obtained |3£] 



ruq > 318 GeV. 



(3.24) 



4. Results 



We now present the results of our study in terms of high relative posterior probability 
regions for CMSSM parameters and superpartner masses (section |4.1| ) and the implications 
for other observables (section ^^ ). In section 4^ we compare those results with the mean 



quality of fit statistics, while in section 4.4 the prospects for direct detection of DM are 
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presented. Correlations among observables are depicted in section 4^, while section 4^ 
is concerned with the influence that the measurement of the anomalous magnetic moment 
has on our conclusions. 

Our statistical inferences are drawn from multiple MC chains, which contain a total 
of about 3 X 10^ to 4 x 10^ samples. For more details about our numerical implementation 
of the MCMC algorithm, see appendix 

4.1 High probability regions for parameters and superpartners masses 

In the six panels of Fig. ^, we show the 2-dimensional posterior relative probability den- 
sity functions p{9i,9j\d), where {6i,9j) = (mi/2,mo), (tan/3,mo), {Ao,mo), {Aq, 1711/2), 
(tan/3, mx/2) ^-nd (tan/3, Aq)) for the "4 TeV range" case which we will treat as our de- 
fault choice. In each panel all other basis parameters have been marginalized over. Redder 
(darker) regions correspond to higher probability density. Inner and outer blue (dark) solid 
contours delimit regions of 68% and 95% of the total probability, respectively. In all the 
2-dimensional plots, the MC samples have been divided into 70 x 70 bins. Jagged contours 
are a result of a finite resolution of the MC chains. 

It is clear that the structure of the parameter space is rather complex. In particular, 
in the left top and bottom panels at mi/2 — 0.2 TeV we can see a narrow high-probability 
funnel induced by the light Higgs boson resonance, which was also observed in the analysis 



of Ref. [17|. The presence of such narrow wedges makes the exploration challenging for the 
MCMC procedure, and much harder for a fixed-grid scan. 

Values of mg < 2 TeV are favored, but larger values are definitely not excluded. In 
particular, the 95% probability region extends up to the upper prior range for mo, a clear 
sign that the data is not powerful enough to constrain this parameter sufficiently (we 
comment further on this issue below). The most probable region in the {mi/2,'m-o) plane 
is centered around the point 

mi/2 ^ 0.7 TeV, mo ~ 0.8 TeV. (4.1) 

The region encompassing 68% of joint probability is roughly bounded by 0.5 TeV < mi/2 
1.5 TeV, because of the efficiency of both the stau coannihilation and/or the pseudoscalar 
resonance. A sharp probability drop above mi/2 ^ 1-5 TeV, almost independent of mo, is 
caused by a combination of the relic abundance and the (Ja^^^"^ constraints. At smaller mo 
the boundary bends because below it the neutralino is not the LSP. Large values of tan /? 
between about 45 and 55 are definitely more favored but smaller values are also allowed, in 
particular for small mo < 0.5 TeV (upper right panel), as a consequence of the light Higgs 
boson resonance. 

The large mo region, starting from the upper part of the 68% probability contour, 
corresponds to a wide range of possible positions of the FP region. For each fixed choice 
of parameters, the FP region consistent with the CDM abundance is actually very narrow, 
but its position varies widely along mo when we marginalize over all the other parameters. 

We observe that Aq is largely uncorrelated with other variables, and its pdf presents 
a strong peak around ^0 — 0.8 TeV. The high probability contours for ^0 are well within 
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tanp tan(3 

Relative probability density 
0.2 ^ ^^^^^^^^^^^1 

Figure 2: The 2-dimensional probability densities in the planes spanned by the CMSSM pa- 
rameters: mo, Ao and tan/? for the "4 TeV range" analysis (see Table |^). The pdf's are 
normalized to unity at their peak. The inner (outer) blue solid contours delimit regions encom- 
passing 68% and 95% of the total probability, respectively. All other parameters in each plane have 
been marginalized over. 
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Figure 3: The 2-dimensional relative probability densities as in Fig. g but for the "2 TeV range" 
analysis and for the planes {mi/2,'mo) and (tan /3, mo) only. Imposing a prior range ttiq < 2 TeV 
modifies substantially the inferences about the high probability regions for mo. 



the prior range (|^o| < 7 TeV), which indicates that this constraint is robust with respect 
to changes to the prior. 

In order to examine the sensitivity of these results to the assumed ranges of CMSSM 
parameters, i.e., the prior used, in Fig. |3| we plot the 2-dimensional pdf's p(W'i/2; ^o|c^) 
and p(tan/?, mo|d) for the "2 TeV range". This is similar to the prior used by the authors 
of Ref. and for this case we find a fairly good general agreement with their results. 
However, several differences in the treatment of uncertainties, in the data employed and in 
the accuracy of the theoretical predictions add up to appreciable differences in the details 
of the results.^ 

It is important to stress that by imposing an upper prior range niQ < 2 TeV one cuts 
away a large region of parameter space which is not excluded by the data (the FP region 
at large mo, compare with the corresponding panels in Fig. ^). Therefore, inferences on 
high probability regions for large mo are different for the "2 TeV range" and the "4 TeV 
range" cases. This means that current data is not informative enough to strongly constrain 
the value of mo independently of prior information, i.e. the prior range one chooses to use. 
The main reason behind this is the presence of the FP region which so far has not been 
investigated by MCMC techniques. On the other hand, results for the other CMSSM 
parameters mi/25tan/3 and Aq do not vary appreciably if one changes the prior ranges, as 
we now discuss. 

In Fig. ^ we show the posterior 1-dimensional pdf p{6i\d) for each of the CMSSM 
parameters, with all the other basis parameters marginalized over. The 1-dimensional 
pdf's contain the complete statistical information about each of the CMSSM variables, 

^For instance, in the analysis of Ref. mo as small as 1 TeV is allowed in the vicinity of the resonance 
(compared to our mo > 1.4 TeV) as a result of employing a less restrictive chargino mass bound of only 
67.7 GeV. 
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Figure 4: The 1-dimensional relative probability densities p{6i\d) for each of the CMSSM param- 
eter, 9i = mo, TOi/2, Aq and tan/3. All other parameters have been marginalized over. The two 
curves compare the results for two different prior ranges (see Table |l|). 



fully accounting for all sources of uncertainty included in the analysis. We plot results 
for both the "2 TeV range" and the "4 TeV range" to facilitate the comparison between 
the two. In the upper left panel, we note that the pdf's for mo < 2 TeV are in excellent 
agreement for both ranges, but above that value the posterior pdf for the "2 TeV range" 
is sharply cut by the prior. On the contrary, the pdf for the "4 TeV range" extends all 
the way to 4 TeV. Since the posterior does not drop to zero before reaching the prior 
range of 4 TeV, it is likely that by allowing even larger values of tjiq one would find non- 
negligible probability densities even there. This is caused by a high sensitivity of the 
position of the FP region along the rriQ axis to the Yukawa couplings ht and h),. The 
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Parameter 


z lev range 
68% region 95% region 


4 lev range 
68% region 95% region 


niQ (TeV) 
mi/2 (TeV) 
Ao (TeV) 
tan/3 


< 1.36 < 1.88 
(0.54,1.20) (0.26,1.48) 
(-0.18,2.01) (-1.44,3.40) 
(44.1,54.9) (16.0,58.0) 


< 2.10 < 3.66 TeV 
(0.52,1.26) (0.20,1.61) 
(-0.34,2.41) (-1.95,4.31) 
(38.5,54.6) (13.6,57.8) 



Table 5: CMSSM parameter ranges corresponding to 68% and 95% of posterior probability (all 
other parameters marginalized over) for the two different prior choices, the "2 TeV range" and the 
"4 TeV range" . 



effect is partially accommodated by treating the top and bottom masses in the nuisance 
parameters. However, as explained in the discussion following Eq. (3.10), the computation 
of /i and as a result of ^x^'^ highly uncertain in the FP region. At present this makes 
it difficult to make a more definitive statement about the region rriQ > 2 TeV other than 
that present data does not strongly constrain it. In fact, the upper 95% probability region 
for the "4 TeV range" extends to mo < 3.66 TeV. For smaller mo, the bulk of the pdf lies 
around 0.5 TeV < mo < 1-5 TeV. 

Constraints on the other CMSSM parameters are largely independent of the adopted 
prior ranges. The bulk of the pdf for mx/2 lies around mi/2 ~ 0.75 TeV, with the 68% 
region within 0.52 TeV < mi/2 < 1-26 TeV, with again a narrow peak due to the light Higgs 
resonance at smaller values. For the "4 TeV range" this narrow peak is more pronounced, 
because in this case one integrates over a larger range for mo, compare with the upper left 
panel of Fig. |2[ In the lower left panel of Fig. ^ the peak around Aq ~ 0.8 TeV is again 
clearly visible and almost independent of the prior used. We notice that the 68% region is 
bounded by —0.34 TeV < Aq < 2.41 TeV and thus Aq = lies close to its boundary. This 
is an interesting result (which can also be seen in Ref. ||l^) in light of the fact that most 
of fixed-grid scans (with a few exceptions |21, Q) have assumed Aq = 0. Finally, in the 
last panel, the 1-dimensional pdf for tan/3 shows a preference for large values, with the 
68% region given by 38.5 < tan/3 < 54.6. Regions containing 68% and 95% of posterior 
probability for the CMSSM parameters are summarized in Table |5|. 

In Fig. H we show 1-dimensional pdf's for several superpartners masses, and the cor- 
responding 68% and 95% probability regions are given in Table ^. Note that the 95% CL 
experimental bounds on the superpartner masses have been included in the likelihood (and 
smeared out by corresponding theoretical uncertainties), as explained in section The 
masses of the gluino g, the lightest chargino Xi^ a-nd the LSP neutralino Xi which are 
proportional mainly to mi/2) are basically the same for both prior ranges, in agreement 
with the result for mi/2 displayed in Fig. ^. In contrast, the pdf's for the masses of the 
sfermions exhibit a sharp cutoff in the "2 TeV range" case, as a consequence of a basic 
mass relation m^r ~ nin + c ? m? ,„ . Therefore the prior cut on mo for "2 TeV range" 

fL,R ^ JL,R 1/2 1- U b 

impacts on the posterior probability distribution for the sfermions as well. Nevertheless, 
for the squarks the relative probability peaks below some 2 TeV which is generally well 
within the LHC reach. As a result, the integrated probability for m^^ < 2.5 TeV is 85%. 
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buper- 


I lev 


range 


4 lev 


^ 77 

range 


partner 


OO70 


95% 


68% 


9570 




Xi 


(0.22, 0.52) 


(0.10, 0.64) 


(0.22, 0.55) 


(0.07, 0.70) 


± 
Xi 


(0.42, 0.98) 


(0.16, 1.20) 


(0.36, 1.00) 


(0.11, 1.25) 


9 


(1.27, 2.64) 


(0.70, 3.19) 


(1.25, 2.80) 


(0.54, 3.51) 


en 


(U.DD, i.by ) 


(U.oU, i.y r ) 


(U.oU, z.yo ) 


(U.oo, o.oo ) 


V 


(0.68, 1.49) 


(0.42,1.76) 


(0.79,2.63) 


(0.46,3.57) 


h 


(0.36,1.03) 


(0.23,1.41) 


(0.42,2.12) 


(0.25,3.31) 




(1.47,2.61) 


(1.10,3.11) 


(1.60,3.50) 


(1.18,4.49) 


ii 


(1.09,2.04) 


(0.82,2.48) 


(1.17,2.44) 


(0.87,3.22) 


hi 


(1.23,2.26) 


(0.98,2.74) 


(1.33,2.79) 


(1.03,3.66) 



Table 6: Selected superpartner mass ranges (in TcV) containing 68% and 95% of posterior prob- 
ability (all other parameters marginalized) for the two different prior choices, the "2 TeV range" 
and the "4 TeV range" . 



For comparison, we find Trig < 2.7 TeV with 78% probability and m ± < 0.8 TeV with 65% 
probability. We will come come back to the issue of the posterior probability distribution 
for the superpartners in subsection ET 



4.2 High probability regions for other observables 

Our MCMC approach allows us to investigate the joint posterior probability distribution 
between CMSSM parameters and the various observables, as explained in section |2.1| . In 
Figs. I - I we plot the joint pdf for Q^h'^, 6af^^^ , BR{B ^^7), BR{Bs 
respectively, and mo, ^1/2 and tan/3 (all other parameters in each panel are marginalized 
over). All plots correspond to the "4 TeV range". (See also column one in Table ^.) We 
stress that the joint pdf is obtained by taking into account observational constraints from 
all the derived variables, including the one plotted along the vertical axis. 

In Fig. |6|, the pdf peaks around mo ~ 1 TeV and ^^h"^ ~ 0.12, with all values of mg 
up to 4 TeV compatible with the observed cosmological DM abundance. This is another 
demonstration that the narrow "WMAP strips" - which appear when including only two 
parameters of the CMSSM at a time (see, e.g., |97, |l5[) - actually widen considerably to 



cover a large region of parameter space when all the variables are taken simultaneously 
into account. Not surprisingly therefore, we do not find a strong correlation between r^^^/i^ 
and the CMSSM parameters. 

In the planes spanned by 5a^^^^ and the CMSSM parameters (Fig. we observe a 
strong anti-correlation between (^a^^^"^ and both mg and mj^/2- This is a simply a reflection 
of the decreasing CMSSM contributions to (5a^^^^ with increasing superpartner masses. In 
general, the posterior distribution for Ja^^^"^ is quite skewed with respect to the likelihood, 
which represents the experimental measurement. The posterior pdf tends to prefer values 
of (5a^^^^ close to zero. We comment further on this below. 

For BR{B Xg'y), Fig. ^ shows a positive correlation with the masses, with the SM 
value being reached in the asymptotic regime of large mi/2 and mo- We also see that the 
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Figure 5: As in Fig. ^ but for the masses of several representative superpartners. 

peak in the pdf is centered slightly below the experimental central value given in Eq. (p.^). 

In Fig. |9|we show the 2-dimensional pdf for BR{Bs — > /U^/i~) and mi/2, rnQ and tan (3. 
As expected, the SUSY contribution decreases with increasing superpartner masses and 
rapidly increases proportionally to tan^ f3. Most of the high relative probability density lies 
close to the SM prediction given in of Eq. dll) 11]. The 1-dimensional (2 tails) regions 
encompassing 68% and 95% of probability are 

3.5 X 10"^ < BR{Bs < 1.68 x lO^^ (68% region), 

3.3 X 10"^ < BR{Bs /u+/i") < 7.50 x 10"^ (95% region). ^^'"^^ 

The current CDF and D0 limits are thus only approaching the 95% probability region. The 
68% (95%) region extends just below (well above) the Tevatron reach of about 2 x 10~^. 
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Figure 6: The 2-dim relative probability density p{il^h'^,6i\d), where 9i — tuq, mi/2 and tan/3. 
Note that the measured value of ^x^^' ^1- ( 3-10 ), has been included in computing the relative 
probability density. 




Note also that a positive measurement of the BR at the Tevatron would imply mo < 
1.5 TeV and would thus strongly disfavor the FP region. 

When combining many different constraints, it is possible that some combinations of 
them might be in conflict with each other. This has been mentioned above for 5a^^^^, 
where we remarked that the posterior pdf tends to prefer a value close to zero. This is 
further highlighted in Fig. IC for a few representative cases. For each observable, we plot the 
1-dimensional marginalized posterior pdf for the two prior ranges, along with the Gaussian 
likelihood function used in the analysis. If all the observations agreed with each other and 
in the absence of strong correlation among variables, we would expect the posterior pdf 
and the likelihood to overlap. This is the case for example for the bottom mass mb(mb)*^'^ 
and for ^^h"^ (upper panels). In the latter case, the slightly skewed shape of the posterior 
is due to our treatment of the theoretical uncertainty, which is larger for larger values of 
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Figure 8: As in Fig. |, but for BR{B Xs^f). The measured value of BR{B X^j), Eq. i^), 
has been included. 
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Figure 9: As in Fig.||, but for BR{Bs /i+M")- The upper limit oiBR{Bs M+A^") < 1.5x 10"'^ 
has been included. 



n^h'^ . However, the posterior pdf for sin'^ 9^^ and 5a^^^^ show a tension with the likehhood 
representing the experimental result. We see a "pull" in the posterior towards values of 
s\v? lower than the measured mean (about one standard deviation discrepancy). We 
notice a similar situation for Mw- This tension is even more pronounced for the SUSY 
contribution to the anomalous magnetic moment of the muon 5a^^^^, whose posterior 
pdf peaks at values close to zero, in contrast to the experimental measurement. This 
is a sign of a tension between the various constraints used, with the other measurements 
pulling the posterior pdf for Ja^^^"^ towards the SM value. This motivates us to investigate 
the dependence of our results on the inclusion of the Ja^^^'^ measurement, which will be 
presented in section |4.6| . 

4.3 Mean quality of fit 

In Bayesian statistics, the posterior pdf represents our state of knowledge about the pa- 
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Figure 10: An illustration of the tensions between different observables. While for mb{mb)^^^ 
and i}^h^ the likelihood and the posterior pdf agree very well (top panels), for sin 9^^^ and (SaSUSY 
(bottom panel) they exhibit a substantial discrepancy. The latter case is a clear sign that the other 
constraints are strongly pulling away the posterior pdf from the value preferred by the anomalous 
magnetic moment measurement. 



rameters in the problem after we have seen the data and given our choice of priors, as we 



have explained in some detail in section 2.1. It is important to remark that regions of high 
posterior probability do not necessarily correspond to the best fitting points. The quality 
of fit is defined in terms of an effective obtained from the likelihood as 

x\9) = -2lnp{d\m)- (4.3) 
We can easily imagine a situation where a tiny multi-dimensional region in parameter space 
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- let us call it a "spike" - shows an excellent quality of fit. At the same time, there might 
be another region where the quality of fit is slightly worse, but whose volume in parameter 
space is much larger. The much larger volume of the latter region gives it a higher weight, 
since it is more generic and hence less fine-tuned. In such a case, the posterior pdf would 
show a strong peak in the large region, while for the spike it would be suppressed due to 
its smallness. At the same time, the quality of fit statistics would show a stronger peak in 
the spike. We emphasize that this is not a feature of the MCMC exploration of parameter 
space, but rather a characteristics built into the Bayesian formalism. As a consequence, 
our inferences automatically give more weight to regions of parameter space showing less 
fine-tuning. 

In the above scenario, an analysis performed using a fixed-grid scan and a quality of 
fit statistics would reach potentially very different conclusions.^ It can be shown, however, 
that both methods lead to the same conclusions if the data is informative enough, i.e., if 
their constraining power is sufficient. Conversely, a discrepancy between the posterior pdf 
and the quality of fit statistics is a useful indicator that the above mechanism is at work. 

In order to carry out such a comparison, we compute the mean quality of fit in two 
dimensions. This is obtained from the posterior pdf by adopting the effective defined in 



Eq. ( [4.3| ) as a function of the parameters f{9), as explained below Eq. (2.7), and plotting 
it in the dimensions of interest in parameter space. In Fig. |ll| we plot the mean quality 
of fit in the same six panels for which we presented the posterior pdf in Fig. ^. The blue 
(solid) contours are the same as in Fig. ^ and are displayed to facilitate the comparison 
between the two quantities. 

In some panels, the best fitting points (represented by dark red regions) are much 
more localized then the high-probability pdf. For example in the {mi/2,'mo) plane we find 
the best quality of fit in the region mo,mi/2 ^ 0.5 TeV. In the (tan /3, mo) plane, on the 
other hand, we observe that good fitting point exist for almost all values of tan /3, down 
to tan/3 ~ 15. Comparing with the corresponding panel in Fig. |2|, we conclude that the 
strong preference for a large tan /? shown by the posterior pdf does not imply that all the 
best fitting points lie in that region of parameter space. This issue can only be resolved 
once better data becomes available. 

4.4 Direct detection of DM 

Predictions for are usually determined as a function of the CMSSM parameters by 
rigidly enforcing relevant constraints, e.g., 1 or 2a ranges of ri^/i^ or BR{B — > ^^7), 
etc. In our analysis, we present a posterior pdf which simultaneously accounts for all the 
constraints and sources of uncertainties. 

In Fig. 12 we plot the posterior pdf for the spin-independent elastic cross section 0"^^ 



and the CMSSM parameters mo, ^-1/2 and tan/3. In the left panel one can see three 
well-separated high probability regions. One is centered at 0.8 TeV < mg < 3 TeV and 
(jp^ ~ 10"^*^ pb, although values almost as large as 10~^ pb are also allowed. It comes from 



®A direct comparison between a Bayesian and a fix-grid (frequentist) analysis would be diflicult, since 
the latter cannot easily handle nuisance parameters. 
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Figure 11: The quality of fit in the planes spanned by the CMSSM parameters: r7ii/2, tiq, and 
tan P for the "4 TeV range" scan. This figure should be compared with Fig. ^ where the relative 
2-dim joint relative probabilities are plotted for the same pairs of variables. 



the bulk, the stau coannihilation and the pseudoscalar resonance regions. To the left, and 
almost adjacent to it, there is a fairly narrow vertical band resulting from a light Higgs 
resonance mentioned above. The last main region is at large mo > 2 TeV and almost 
constant a^^ ~ 1.6 x 10~^ pb which results from the FP region. Comparing with the right 
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panel in Fig. 12 we can see that dp from the FP region and the light Higgs resonance 
are fairly independent of tan /3 while the bulk, coannihilation and resonance region requires 
large tan f3. Finally, the middle panel is added for completeness and it closely resembles the 
top panel of Fig. ^ (because ~ OAmi/2) where we present the 2-dimensional pdf for 
CTp^ and m^. For comparison, we also show current CDMS-II, Edelweiss-I and UKDMC 
ZEPLIN-I 90% CL limits, but we stress that this constraint has not been used in the 
likelihood. 

The features discussed above are visible in Fig. |13| (top). Firstly, the biggest, banana- 
shaped region of high probability (68% regions delimited by the internal solid, blue curve) 
shows a well-defined anti-correlation between cr^^ and m^. It comes from the bulk and 
stau coannihilation region (larger 0"^^) and from the ^-resonance (smaller (Jp^) at large 
tan/3. This region covers roughly the range 10^^*^ < < 10^^ pb and 0.2 TeV < < 
0.7 TeV. In both cases the dominant contribution to 0"^^ typically comes from a heavy 
Higgs exchange. 

Secondly, at small ^0.1 TeV we can again see a small vertical band of fairly low 
relative probability density (< 0.2) at small 0"^^. This region is allowed by the light Higgs 
resonance contribution to reducing f^^^^ small mi/2- This region essentially never 
features in usual fixed-grid scans, which do not smear out experimental bounds. This 
region would also disappear with a fair improvement in the lower bound on m/j. 

Thirdly, we can see a well pronounced region of high relative probability at fairly con- 
stant (Tp^ ~ 1.6 X 10~^pb for iTLy^ < 0.42 TeV which at low partly overlaps with the 
previous region. At 95% this region extends up to < 0.72 TeV for fairly constant cip^. 
This "high" a^^ band is a result of the FP region, basically independently of tan /?, as 
discussed above. This result has to be interpreted carefully, since there are large uncer- 
tainties associated with FP region, in particular with its location in the {1711/2,^0) plane 
mentioned earlier. Hopefully, associated uncertainties in 0"^^ are going to be much smaller 
since it depends on low energy quantities like Higgs masses and the fi parameter. Despite 
those outstanding questions, we believe that it is probably safe to expect that the FP will 
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Figure 13: Top panel: the 2-dimensional relative probability density p{m^,ap^\d), with the 
contours containing 68% and 95% probability also marked. Bottom panel: the mean quality of fit 
(likelihood) with the same 68% and 95% probability contours as in the top panel to facilitate the 
comparison. Current 90% experimental upper limits are also shown. 



be the first to be probed by DM search experiments. 

Finally, after marginalizing over all other parameters, we obtain the following 1- 



-33- 




Relative probability density 
0.1 ^^0.2 0.3 0.4 O ^j^^^^^^^^^^^^^^l 



Figure 14: The 2-dimensional relative probability density p{^i,£,j\d) for the pairs of observables 
marked along horizontal and vertical axes. 



dimensional regions encompassing 68% and 95% of the total probability: 

1.0 X 10"^° pb < cr^^ < 1.0 X 10"^ pb (68% region), 

0.5 X 10"^° pb < a^^ < 3.2 x 10"^ pb (95% region). ^^'^'^ 

Currently running experiments (most notably CDMS-II but also Edelweiss-II and 
ZEPLIN-II) should be able to reach down to a few x 10~^ pb, on the edge of exploring the 
FP region. A future generation of "one-tonne" detectors are expected to reach down to 
(Tp^ > lO^^'^ pb, thus exploring almost the whole 68% region and much of the 95% interval 
as well. 



As discussed in section 4.3, the posterior pdf may be rather different from the mean 



quality of fit. In the bottom panel of Fig. 13 we plot the quality of fit (defined in Eq. (|4.3|)) 



for (Tp^ and m^. The best fit points are found in the region 0.1 TeV < < 0.2 TeV and 
1 X lO^^'^ pb < ^ 3 X 10^^ pb, but other good fitting points (quality of fit about 0.4, 
dark green regions) lie right near the top of the high probability region, at rather large 



Up 



4.5 Correlations among observables 

We now proceed to examine various correlations among the observables discussed above. 
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In the upper left panel of Fig. 14 we show the 2-dimensional pdf for BR{Bs ^ ) and 

can clearly see a rather strong correlation between those two variables |^8| 
but, as pointed out above, the most probable ranges of both variables are on a low side. A 
positive measurement of BR{Bs fi^fi^) at the Tevatron (above some 2 x 10^^) would 
lead to a very strong tension with the current result for Ja^^^'^. 

In the upper right panel of Fig. ^ we show the 2-dimensional pdf for BR(Bs — > jj,^ 
and cTp^. We can see two high relative probability regions here showing an interesting 
pattern. The one at smaller values of both variables comes from the coannihilation and 
Higgs resonance regions and has been pointed out in |99|. In addition, the FP region 
allows for a second "island" where a^^ is not far below the current CDMS limit while 
BR{Bs — > /i^/i^) is very small, far below the reach of the Tevatron. Thus, because of the 
FP region, a signal in measuring 0"^^ in current generation of DM search detectors would 
not necessarily imply a high probability of measuring BR{Bs /U^^^) at the Tevatron, 
(or vice versa, contrary to the claim of p9t|). 



In the two bottom panels of Fig. we show high probability regions for BR{Bs 
fi'^fi") and BR{B — > Xg^), and for BR{B —>■ Xg^) and 0"^^. In the first case, the high 
relative probability region for BR{Bs fJ-^fJ- ) lies at small values, just above the SM 
prediction and corresponds to the BR{B — > Xg'j) around the central value, as noted 
above. It may eventually be possible to verify this correlation experimentally. Finally 
the variables BR{B Xgj) and a^^ show two separate peaks in the relative probability 
density. Again, the peak at smaller cXp^ comes from the coannihilation and Higgs resonance 
regions while the other one from the FP region. In principle, a measurement of above 
some 10~^ pb would point towards a range of BR(B — > Xg'j) above the current central 
value. Unfortunately, until a full NLO SUSY contribution is available, error bars in the 
latter quantity will remain at the level of some 10%, which would make it difficult to 
confirm the FP origin of the 0"^^ measurement. 

In the four panels of Fig. |l^ we plot in the {mi/2,'m'o) plane values of rrih (upper left), 
^x^'^ (upper right), (5a^^^^ (lower left) and cr^^ (lower right). The points have been drawn 
uniformly from our MC chains. To highlight the values of interest for the observables, 
the range of the color scales has been reduced, and points with values above (below) the 
scale have been plotted in red (blue). One can see that nih increases with increasing mx/2 
or mo, as expected. In the region where mi/2 ^ 1 TeV and/or rriQ < 2 TeV (roughly 
the reach of the LHC) one predominantly finds ruh ^117 GeV (although larger values 
are not excluded), which is encouraging for Higgs searches |10C1[ |. In all the panels, at 
mi/2 — 0.2 TeV, there is a vertical favored region due to a narrow light Higgs resonance 
contribution to It is interesting that in the rest of the (m-1/21 "t-o) plane one finds 

that the WIMP relic density (upper right window) can take any value within about the 
2a range of ( 3.12 ). In other words, even though for each particular choice of the CMSSM 
and nuisance SM parameters there are only a few narrow regions consistent with the DM 
constraint (3.10), by performing the appropriate marginalization over all other parameters 
it appears to be fairly easy to satisfy the WMAP value of the DM abundance at almost 
any point in the (?t^i/2)™'o) plane. This feature can also be seen in Fig. |^. 
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Figure 15: Values of rrih (upper left), Q^h"^ (upper right), a^^ (lower left) and (5a^^^^ (lower 
right) in the (mi/2,'Tio) plane for samples uniformly selected from our MC chains. While ct^^ and 
(5a^j^^^ show a clear correlation with the values of (mi/2,nT-o), the CDM abundance ^1x1^-^ can take 
any value within the 2a WMAP range in the plane. Narrow "WMAP strips" in the (mi/2,'Tio) 
plane for fixed values of tan/3 and Aq thus disappear when these parameters are allowed to vary 
and all other parameters are correctly marginalized over. 



The values of a'^^ (bottom left panel of Fig. 15) cover a large range but are generally 
larger for smaller mi/2 and intermediate and large mo (compare with Figs. |l^ and |l^ 
(top)), with the exception of the light Higgs boson resonance region at small mi/2- This 
means that DM direct detection searches will in general explore the large mo region of FP 
first ~ an interesting complementarity with collider searches. 

On the other hand, the claimed discrepancy between experiment and the SM value ( |3.7| ) 
of the anomalous magnetic moment, when taken at face value, clearly points towards a 
different region of mi/2 0.8 TeV (at small mo) and mo < 1.5 TeV (at small mi/2), as can 
be seen in the bottom right panel of Fig. |l^. Thus, similarly to the case of BR{Bs — > 
described above, a positive measurement of a^^ in currently running DM detectors (above 
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Figure 16: The 1-dimensional relative posterior pdf 's as in Fig. ^. The black solid (blue dotted) 
line is for the analysis including (excluding) the measurement of the anomalous magnetic moment 
of the muon, Eq. (3.7). Both cases assume the "4 TeV range" priors. 



some 2 X 10 ® pb) would lead to a strong tension with the current result for (5a^^^^. 
4.6 Sensitivity to {g — 2)^^ 

As we have shown above, the 2.7a deviation of the anomalous magnetic moment of the 
muon from the SM prediction appears to be in tension with some of the other observables. 
To investigate to what extent our statistical inferences on superpartner masses rely on the 
inclusion of (5a^^^^, in this section we present the posterior pdf's obtained by dropping it 
from the likelihood. To be as general as possible, we have adopted the "4 TeV range" of 
the priors. 
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Figure 17: As in Fig. ||. The black solid (blue dotted) line is for the analysis including (excluding) 



the measurement of the anomalous magnetic moment of the muon, Eq. (3.7). Both cases assume 
the "4 TeV range" priors. 



Fig. 16 compares the 1-dimensional marginalized posterior pdf p{Oi\d), where 9i = uiq, 
^0 and tan f3 with and without the inclusion of (5a^^^^ in the likehhood. If we drop 
the anomalous magnetic moment measurement we loose essentially all constraint on uiq, 
whose pdf becomes flat above 1 TeV. The impact on mi/2 is very mild, with a slight shift 
to larger values of the bulk of the pdf. There is also almost no change in the pdf 's for 
and tan /3. Intervals encompassing 68% and 95% of probability for the CMSSM parameters 
are given in Table ^ (compare with Table ^). 

The implications for several representative superpartner masses are presented in Fig. 17 , 
which should be compared with Fig. ^. The corresponding probability intervals are sum- 
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Figure 18: The probability for the superpartner masses to lie below a given mass. All other 
parameters have been marginalized over. 



marized in Table ^. One should remember that the dominant contribution to 5a^^^^ comes 



from the sneutrino-chargino exchange llOlj. Since m-p depends much more strongly on mo 
than on mi/2, while '^^i is more dependent on mi/2, both soft parameters are affected 
but rriQ more strongly. 

An experimental implication for the LHC is rather obvious. If the 5a^^^^ anomaly 
is not confirmed, the probability of finding sleptons and squarks will be reduced. This 
can be see in Fig. |l^ where we plot the total probability as a function of mass for several 
superpartners. We show three cases: the "2 TeV range" (dotted red), as well as our 
default "4 TeV range" case with (black solid) and without (dotted blue) the {g — 2)^ 
constraint ( |3.7D . It is clear that the total probability that the mass of squarks or sleptons 
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Parameter 


4 lev range , no [g - 2)^ 
68% region 95% region 


mo (TeV) 
mi/2 (TeV) 
^0 (TeV) 
tan/3 


< 2.87 < 3.85 
(0.56,1.38) (0.14,1.73) 
(-0.60,2.89) (-2.46,4.59) 
(28.2,53.8) (10.2,57.0) 



Table 7: CMSSM parameter ranges corresponding to 68% and 95% of posterior probability (with 
all other parameters marginalized over) for the "4 TeV range" but without the constraint from 
{g — 2)fj,. These ranges should be compared with Table ^ 



Super- 


"4 TeV ran^ 


5e", no (5-2);, 


partner 


68% 


95% 




(0.23,0.60) 


(0.57,0.76) 


xt 


(0.33,1.10) 


(0.11,1.31) 


9 


(1.36,3.05) 


(0.42,3.79) 




(1.12,3.48) 


(0.40,3.95) 


V 


(1.11,3.18) 


(0.56,3.75) 


n 


(0.61,2.71) 


(0.32,3.58) 




(2.08,4.00) 


(1.43,4.75) 


ii 


(1.47,2.76) 


(1.00,3.48) 


bi 


(1.70,3.19) 


(1.21,3.95) 



Table 8: Selected superpartner mass ranges (in TeV) containing 68% and 95% of posterior prob- 
ability (with all other parameters marginalized) for the "4 TeV range" but without the constraint 
from (5 — 2)^. These intervals should be compared with Table ^. 



Observable 


"4 TeV range" , with {g - 2)^ 
68% 95% 


"4 TeV range", w/o {g - 2)^ 
68% 95% 


5aSUSY X lOio 

BR{Bs /i+Ai") X 10^ 
BR{B Xs-f) X 10^ 
cj^^(pb X 10^°) 


(0.107,0.138) (0.946,0.157) 
(1.9,9.9) (0.8,17.1) 
(3.5,16.8) (3.3,75.0) 
(2.93,3.44) (2.62,3.61) 
(1,100) (0.5,320) 


(0.107,0.138) (0.949,0.158) 
N/A N/A 
(3.5,8.5) (3.3,41.0) 
(3.08,3.49) (2.77,3.62) 
(0.8,117) (0.3,208) 



Table 9: Intervals encompassing 68% and 95% of posterior probability for selected observables 
(with all other parameters marginalized). The two columns compare results for the "4 TeV range" 
with (left column) and without (right column) inclusion of the constraint from {g — 2)^. 



lies below a certain value depends rather strongly on the choice priors and/or the {g — 2)^ 
constraint. On the other hand, this is not the case for the gauginos. By comparing Tables ^ 
and ^ it would seem that the upper bound of the 95% interval does not change much for 
the squark and gluino masses. However, this is probably a manifestation of the upper cut 
induced around 4 TeV by the prior, and therefore should not be taken as a robust result 
of the inference, analogously to what has been discussed above for niQ. 
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In Table ^ we summarize the effect of the {g — 2)^ constraint on some of our other 
observables. While the relic abundance is unaffected, BR{Bs — > iJi^ IjT) and a^^ get shifted 
to lower values as a consequence of a less constrained mo when the [g — 2)^ anomaly is 
discarded. 

5. Summary and conclusions 

In this work we performed a detailed investigation of the CMSSM parameter space using 
a MCMC method and analyzed our results in terms of Bayesian statistics. The power 
and flexibility of the approach allowed us to probe many previously unexplored ranges 
of parameters. Furthermore, we were able to improve on previous analyses in several 
important aspects. We fully incorporated the effects of remaining uncertainties in relevant 
SM parameters (including aem(-^z)^^'^)) which are usually fixed to their central values. 
We incorporated often neglected theoretical uncertainties in computing mass spectra and 
observables. Finally, we improved upon the usual practice of applying sharp experimental 
limits and (typically) la uncertainties by smearing them out. 

By carefully constructing the likelihood function we constrained the CMSSM param- 
eters by comparing several collider and astrophysical observables (except for (t^^) with 
the available data. For all these variables we computed the posterior probability density 
functions ~ our main tool in analyzing our findings - in terms of which we delimited the fa- 
vored regions of the CMSSM parameters, and further derived the ranges of the observables 
themselves that are favored by a combination of currently available data. We emphasized 
the difference between the posterior probability density and the quality of fit statistics. 
The latter is formulated in terms of an average effective x^j and is akin to what is used in 
fixed-grid scans. For example, we concluded that the strong preference for a large tan (3 
shown by the posterior pdf does not imply that all the best fitting points lie in that region 
of parameter space. This issue can only be resolved by acquiring better data. 

We explored in detail the robustness and sensitivity of our results to the assumed 
initial ranges of CMSSM and SM parameters (priors). To this end we compared our 
findings for the default "4 TeV range" prior (which include the somewhat uncertain FP 
region) extending above the LHC reach, with the more restrictive "2 TeV range" prior. 
We emphasized that much care must be exercised in interpreting our inferences whenever 
boundaries of high probability regions lie close to the prior ranges. This applies mainly to 
mo, and to the superpartner masses that primarily depend on it, while all other variables 
appear robust to a change in the range of the priors. 

We furthermore examined various correlations among the relevant observables. Some 
have been pointed out before, e.g., between BR{Bs — > //+//-) and fcSUSY 

or (the last 

one showing a new feature due to the presence of the FP region) . A more subtle correlation 
between BR{B — > Xg'y) and BR[Bs — > fi^ emerged, which may eventually be tested 
experimentally. We note that at present none of the observables appear to be in conflict 
with observations or with each other, with the possible exception of {g — 2)^. In particular, 
the cosmological constraint on Jl^/i^ appears less severe than what has been previously 
thought. 
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Our findings in tlie CMSSM strongly support the idea of low energy SUSY, although of 
course they do not disfavor other possibilities, like split SUSY |102| ]. Quantitatively, at 68% 
probability, we have found 0.52 TeV < mi/2 < 1.26 TeV, mo < 2.10 TeV, -0.34 TeV < 
Aq < 2.41 TeV and 38.5 < tan /? < 54.6. The corresponding ranges of superpartner masses 
are given in column two of Table ^. A significant fraction of the 68% probability ranges of 
superpartner masses falls within the LHC reach, and typically outside the Tevatron reach. 
The same applies to BR{Bs jj,'^ for which most favored CMSSM values are not far 
above the SM prediction. On the other hand, a positive measurement of BR{Bs — > 
at the Tevatron would strongly disfavor large mg, including the focus point region. The 
WIMP DM direct detection elastic scattering cross section a^^ shows a wide spread of 
values (below today's limits) at around 10~^^^ pb and a strong anti-correlation with m^. 
In addition, a relatively large a^^ ~ 1.6 x 10~^ pb, fairly independent of m^, appears to 
be a feature of the FP region (despite large theoretical uncertainties) and will probably be 
the first to be tested in direct detection experiments. 

The {g — 2)^ anomaly still remaining the subject of some controversy, we re-examined 
its impact on the CMSSM parameter space. We showed the inference to be substantial 
on mo, and any superpartner masses that primarily depend on it, while much less so with 
the other CMSSM parameters. The chance for the LHC to detect superpartners reduces 
somewhat but still remains strong. 
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A. Markov chain Monte Carlo algorithm 
A.l Sampling 

The purpose of the Markov Chain Monte Carlo (MCMC) algorithm is to construct a 
sequence of points in parameter space (called "a chain") whose density is proportional 
to the posterior pdf of Eq. (^I5|). Once such a chain has been produced, the posterior 
probability for a given region of parameter space (a bin) is obtained by simply counting 
the number of samples within that region. Marginalization over nuisance parameters (see 

^Available from cosmologist . inf o. 
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Eq. (|2.q )) is trivial: the coordinates of the parameters that one is not interested in are 
simply ignored when counting the samples. 

Several algorithms are available to construct Markov chains, which are more or less 
suited to the structure of the parameter space under consideration, see e.g. |103] for an 



introductory review and references therein. We make use of the Metropolis algorithm: from 
a starting point in parameter space rjQ with associated posterior probability poi a candidate 
point 7/c for the next sample is proposed by sampling it from a transition probability 
T{rio,r]c). The candidate sample is then accepted with probability 

Q = min ( — , 1 ) , (A.l) 

\Po J 

where Pc = p{r]c\d). Notice that all steps for which the candidate sample has a better 
probability than the previous one are accepted. If the candidate point is accepted, it 
becomes the new starting point, and another candidate is drawn. Otherwise the chain 
stays at the previous point (which is thus counted twice) and a new attempt is made from 
there. It can be shown that the sequence of samples r/o, ?7i, • • • , . . . constructed this way 
converges to the target distribution p{r]\d) for t — > oo. 

The transition probability T[r]t,r]t+i) must be symmetric in its arguments (so called 
"detailed balance" ) , a sufficient condition which ensures that the Markov chain constructed 
this way is sampling from the target probability distribution. In our case, the transition 
probability is given by the following prescription, which we found strikes a good balance 
between efficiency and robustness of the exploration. At each step, we alternatively update 
the value of the CMSSM parameters 6 or of the nuisance parameters ij:. In general, all of 
the components of either 9 or are updated at each step. The candidate point is proposed 
along the direction w, where the vector w is restricted to either one of the two subspaces 
(CMSSM or nuisance parameters) and it is given by 

w = A-u. (A.2) 

Here, A is a random rotation matrix which is restricted to a given subspace and which is 
renewed every 4 steps. The components of u are initially chosen as a guess of the typical 
spread of the posterior distribution along each direction of parameter space. The results of 
a preliminary MCMC run are then used to estimate the covariance matrix for the posterior 
pdf, whose eigenvectors give directions of approximate degeneracies in the parameter space. 
In the final run, w is built analogously as above, but this time by a random rotation in 
the space spanned by the projections of the eigenvectors of the covariance matrix. This 
procedure aims at aligning the directions of the proposals to degeneracy lines in parameter 
space, thus improving the efficiency of the MCMC walk. 

Along the direction defined by w, the width of the step is chosen by multiplying \w\ 
by a scaling factor s and a factor r drawn from the distribution 

2 1 

p{r) oc — r"~^ exp(— nr^/2) H — exp(— r), (A.3) 

3 3 



with n = 4 and s = 2.4. The first term on the rhs of the proposal distribution ( [A.3| ) tends 
to make the chain move away from r = for n > 1, while the second term increases the 
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probability near the origin. Thus this distribution tends to be robust even in situations 
where the target pdf is strongly non-Gaussian. The above choices for the proposal distri- 
bution and for the parameters n, s are mostly a matter of trial-and-error. Our updating 
procedure follows closely the recommendations of \W4]. Notice that since our choice of the 
step direction and size do not depend on rjt,r]t+i at any moment, the condition of detailed 
balance for the transition probability holds true. 



A. 2 Convergence 

We start = 12 or = 16 chains in randomly chosen points of parameters space (within 
the boundary specified by our prior range), making sure that they lie well apart from each 
other in order to maximise the initial variance. A certain number of samples have to be 
discarded at the beginning of the chain, since the chain requires some time to equilibrate 
around the target distribution. This "burn-in period" is assessed by inspecting the evo- 
lution of the probability as a function of the number of samples. We find that discarding 
10'^ initial samples is more than sufficient. 

The acceptance rate is defined as the percentage of accepted proposals. For the runs 
using the covariance matrix the typical acceptance rate is in the range 2 — 3%. This is 
rather low, compared to optimal cases where the acceptance rate is typically an order of 
magnitude larger. The reason for this is the configuration of the posterior pdf in multi- 
dimensional parameter space, which is strongly non-Gaussian and presents a challenging 
combination of wide regions and narrow wedges of high probability. There is no optimal 
strategy in this case, but we are confident that our chains have appropriately sampled the 
whole parameter space. We performed extensive checks by comparing runs obtained with 
and without the covariance matrix in order to make sure that the efficiency gain did not 
come at the expenses of a reduced sampling ability. 

Mixing and convergence of the chains is assessed with the Gelman & Rubin R- 



statistics | 105 |. This represents the variance of the means divided by the mean of the 
variances between different chains. The usual criterion is that R — 1 < 0.2, but to be 
conservative we require that for our chains -R — 1 < 0.05 for all parameters (this means 
that our convergence criteria are more stringent). 

Typically we run in parallel two sets of = 12 or = 16 chains, until each chain 
within the set has reached 3 x 10^ or 4 x 10^ samples^ (the exact numbers depending on 
the computing power available and on the convergence status). We check that each run 
has converged using the criterion outlined above, and we then perform consistency checks 
between the two runs. The final inferences are obtained after merging the two sets of 
chains together. At this final stage, the Gelman and Rubin criterion is again satisfied by 
the merged set consisting of 24 < A^ < 32 chains, containing a total number of samples in 
the range of 0.7 x 10^ to 1.3 x 10*^. 

* Since each sample is obtained with a typical acceptance rate of 3%, this means that each chain requires 
0(10'') likelihood evaluations, each of which takes about 1-2 sec on our machines. 
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